function [schidx,nsidx] = SelfTrainSelectFcn(sch,schidx)
% [schidx,nsidx] = SelfTrainSelectFcn(sch,schidx)
%
%  Select the next trial index (NSIDX) based on previous responses. Once
%  n consecutive hits are received in a row, then increase stimulus onset
%  delay by x milliseconds.  If n consecutive misses are received in a row,
%  then decrease stimulus onset delay by x milliseconds.
%
%  Use the following parameters in schedule.writeparams including *
%  ('Write' or 'Read/Write' functionality in ScheduleMaker):
% 
%       *FreqList  ... list of frequencies to test

%  Takes the schedule and schidx as input.
%       SCH is the schedule structure
%       SCHIDX is a 1D matrix with the number of rows equal to the number
%           of unique trials; i.e. size(trials,1).
%  Returns an updated version of SCHIDX and NSIDX
%       NSIDX is the next schedule trial index, i.e. the next index
%           selected from the SCH.trials matrix to run.
%       SCHIDX should be updated with the number of runs for each unique
%           trial in SCH.trials.  The returned version should have
%           SCHIDX(NSIDX) incremented.
%
%  This same format for input and output paramters must be met for custom
%  trial select functions (i.e., experiment.trialselectfcn);
%
%  DJS (c) 2012

global G_RP

persistent TR % tracking

% some constants ----------------------------------------------------------
START_DELAY = 0; % milliseconds
MIN_DELAY   = 0;
MAX_DELAY   = 7000;
INCR_DELAY  = 100;
DECR_DELAY  = 100;
INCR_AFTER_N = 5; % increase delay after N consecutive hits
DECR_AFTER_N = 5;
PERF_CRITERION = 0.8; % performance criterion
RESP_WIN     = 2000; % Response window duration
% -------------------------------------------------------------------------


% CATCH_P = 20; % catch trial %

if ~isfield(sch,'trial_idx'), sch.trial_idx = 1; end
trials = sch.trials;
% boxid  = sch.boxid;
RV = sch.response_vals;

% find some paramters
rcidx = find(ismember(sch.readparams,'response_code'))+1;
% ttidx = find(ismember(sch.readparams,'trial_type'))+1;
sdidx = find(ismember(sch.writeparams,'*START_DELAY'));

if ~isempty(sdidx)
    START_DELAY = trials{1,sdidx};
end

% which module?
writemod = sch.writemodule(rcidx);
RP = G_RP(writemod);

% sampling rate
Fs = RP.GetSFreq;

% convert constants milliseconds -> samples
START_DELAY = round(START_DELAY / 1000 * Fs);
MIN_DELAY   = round(MIN_DELAY   / 1000 * Fs);
MAX_DELAY   = round(MAX_DELAY   / 1000 * Fs);
INCR_DELAY  = round(INCR_DELAY  / 1000 * Fs);
DECR_DELAY  = round(DECR_DELAY  / 1000 * Fs);
RESP_WIN    = round(RESP_WIN    / 1000 * Fs);

if START_DELAY == 0, START_DELAY = 1; end % start delay must be atleast 1 sample

% stimulus duration
sdidx = ismember(sch.writeparams,'stim_dur');

% initialize first trial
if isempty(RV)
    TR = 0;
    
    schidx = zeros(size(trials,1),1); % initialize SCHIDX
    nsidx = 1;
    schidx(nsidx) = 1;
    
    stim_dur = sch.trials{nsidx,sdidx};
    
    RP.SetTagVal('trial_type',0);
    RP.SetTagVal('stim_delay',START_DELAY);
    RP.SetTagVal('init_grace',START_DELAY);
    RP.SetTagVal('stim_end_val',START_DELAY+stim_dur);
    RP.SetTagVal('respwin_end_val',START_DELAY+stim_dur+RESP_WIN);
    return
end

% % Test if catch trial
% tt = rand(1) < CATCH_P;
% RP.SetTagVal('trial_type',tt);

% response codes
RC = RV(:,rcidx);

% delay
DELAY = RP.GetTagVal('stim_delay');

% recent consecutive hits
if TR < 1 && length(RC) >= INCR_AFTER_N
    ind = bitget(RC(end-INCR_AFTER_N+1:end),6);
    if sum(ind) / numel(ind) >= PERF_CRITERION
        DELAY = DELAY + INCR_DELAY;
        fprintf('Trial#: %d\tDelay extended to %0.2f seconds\n',length(RC),DELAY/Fs)
        TR = INCR_AFTER_N;
    end

% recent consecutive misses
elseif TR < 1 && length(RC) >= DECR_AFTER_N
    ind = bitget(RC(end-DECR_AFTER_N+1:end),5);
    if all(ind)
        DELAY = DELAY - DECR_DELAY;
        fprintf('Trial#: %d\tDelay decreased to %0.2f seconds\n',length(RC),DELAY/Fs)
        TR = DECR_AFTER_N;
    end
    
else
    TR = TR - 1;
end

% add some randomization to the delay
a = MIN_DELAY; b = DELAY;
DELAY = rand * (b - a) + a;

% keep DELAY within bounds
if DELAY > MAX_DELAY, DELAY = MAX_DELAY; end
if DELAY < MIN_DELAY, DELAY = MIN_DELAY; end

% select random trial
% give priority to least chosen trials
i = min(schidx);
i = find(schidx == i);
r = randperm(length(i));
x = i(r(1));
schidx(x) = schidx(x) + 1;
nsidx = x;
stim_dur = sch.trials{nsidx,sdidx};

% update parameters
RP.SetTagVal('trial_type',0);
RP.SetTagVal('stim_delay',DELAY);
RP.SetTagVal('init_grace',DELAY);
RP.SetTagVal('stim_end_val',DELAY+stim_dur);
RP.SetTagVal('respwin_end_val',DELAY+stim_dur+RESP_WIN);









